Skip to content

Demand Response Optimization for Battery Energy Storage#679

Draft
vijay092 wants to merge 7 commits intoNatLabRockies:developfrom
vijay092:peakload-optimized
Draft

Demand Response Optimization for Battery Energy Storage#679
vijay092 wants to merge 7 commits intoNatLabRockies:developfrom
vijay092:peakload-optimized

Conversation

@vijay092
Copy link
Copy Markdown
Collaborator

@vijay092 vijay092 commented Apr 15, 2026

Demand Response Optimization for Battery Energy Storage (Stage 1)

This PR introduces a Pyomo-based formulation for demand response, which will be implemented in two stages.

As the first stage, this work implements a rolling horizon optimization for battery operations. The battery dispatch logic is based on a pre-defined signal, such as LMP, load, or a combination of both. This is the G&T level dispatch signal for demand response. The next stage will implement the peak load management logic.

Section 1: Type of Contribution

  • Feature Enhancement
    • Framework
    • New Model
    • Updated Model
    • Tools/Utilities
    • Other (please describe):
  • Bug Fix
  • Documentation Update
  • CI Changes
  • Other (please describe):

Section 2: Draft PR Checklist

  • Open draft PR
  • Describe the feature that will be added
  • Fill out TODO list steps
  • Describe requested feedback from reviewers on draft PR
  • Complete Section 7: New Model Checklist (if applicable)

TODO:

  • Add unit tests
  • Add documentation

Type of Reviewer Feedback Requested (on Draft PR)

Structural feedback: Is this in the right place?

Implementation feedback: I used the same style as Gen's pyomo implementation. Appreciate feedback here.

Other feedback:

Section 3: General PR Checklist

  • PR description thoroughly describes the new feature, bug fix, etc.
  • Added tests for new functionality or bug fixes
  • Tests pass (If not, and this is expected, please elaborate in the Section 6: Test Results)
  • Documentation
    • Docstrings are up-to-date
    • Related docs/ files are up-to-date, or added when necessary
    • Documentation has been rebuilt successfully
    • Examples have been updated (if applicable)
  • CHANGELOG.md
    • At least one complete sentence has been provided to describe the changes made in this PR
    • After the above, a hyperlink has been provided to the PR using the following format:
      "A complete thought. [PR XYZ]((https://github.com/NatLabRockies/H2Integrate/pull/XYZ)", where
      XYZ should be replaced with the actual number.

Section 3: Related Issues

Section 4: Impacted Areas of the Software

Section 4.1: New Files

  • Main Implementation: h2integrate/control/control_strategies/storage/plm_optimized_storage_controller.py

  • Usage Example: examples/34_plm_optimized_dispatch

Section 4.2: Modified Files

  • h2integrate/core/supported_models.py
    • Added the new model

Section 5: Additional Supporting Information

Section 6: Test Results, if applicable

Section 7 (Optional): New Model Checklist

  • Model Structure:
    • Follows established naming conventions outlined in docs/developer_guide/coding_guidelines.md
    • Used attrs class to define the Config to load in attributes for the model
      • If applicable: inherit from BaseConfig or CostModelBaseConfig
    • Added: initialize() method, setup() method, compute() method
      • If applicable: inherit from CostModelBaseClass
  • Integration: Model has been properly integrated into H2Integrate
    • Added to supported_models.py
    • If a new commodity_type is added, update create_financial_model in h2integrate_model.py
  • Tests: Unit tests have been added for the new model
    • Pytest-style unit tests
    • Unit tests are in a "test" folder within the folder a new model was added to
    • If applicable add integration tests
  • Example: If applicable, a working example demonstrating the new model has been created
    • Input file comments
    • Run file comments
    • Example has been tested and runs successfully in test_all_examples.py
  • Documentation:
    • Write docstrings using the Google style
    • Model added to the main models list in docs/user_guide/model_overview.md
      • Model documentation page added to the appropriate docs/ section
      • <model_name>.md is added to the _toc.yml

Copy link
Copy Markdown
Collaborator

@elenya-grant elenya-grant left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi Sanjana! Just putting a few high-level comments now - but I will do a deeper dive shortly! I love the thorough testing of PLMOptimizedControllerConfig and PLMOptimizedStorageController! The docstrings also look awesome! Besides the comments and questions that I left, other high-level things I'd like to see included before this is merged in are:

  1. integration tests of the controller with the StoragePerformanceModel - similar to the tests in h2integrate/control/control_strategies/storage/test/test_optimal_controllers.py
  2. a test for the example you added in examples/test/test_all_examples.py (this can be a pretty basic test, but it'll be good to ensure that the example you added is updated and runs without error as H2I is further developed)

Truly awesome work! Please message or call me if you want to discuss any of the comments I left! I will do another review that focuses on the logic within the controller - but wanted to get these initial questions and comments to you early on!

discharge_efficiency: float = field(validator=range_val(0, 1), default=1.0)
n_max_events: int = field(default=10)
n_control_window: int = field(default=24 * 30) # one month of hourly data
signal_threshold_percentile: float = field(default=0.0)
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Small comment - may be nice to add a range_val validator for the signal_threshold_percent (unless it can be over 100)?

Suggested change
signal_threshold_percentile: float = field(default=0.0)
signal_threshold_percentile: float = field(default=0.0, validator=range_val(0,100))

This is a suggestion/thought/nitpick (but not necessary) - sometimes it's nice to use fraction instead of percentages (specifically with range_val validators) because it is easier to catch with validators.

Aka - if a parameter is defined as a percent (between 0 and 100), but a user thinks its defined as a fraction (suppose they want it to be 40%, but set the value to 0.4 because they think its a fraction), then this would not throw an error with a validator. But - if a parameter is defined as a fraction (between 0 and 1) and a user thinks its defined as a percentage (they want it to be 40% so they set the value to 40), then the validator does thrown an error.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Typically, percentiles are not written in terms of fractions (unlike percentages). Have you seen this done before?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also, that sounds good about the range_val. I hadn't added it because np.percentile will catch it anyway.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it would be worth including the validator even if np.percentile will catch it later because the validator will provide a very clear error message and appear readable for users. non-blocking comment.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That sounds good, I added it.


max_charge_rate: float = field()
supervisory_signal: list = field()
peak_window: dict = field()
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a limit/constraint on peak_window specifically with respect to n_control_window? Like what if someone defined n_control_window = 6 and peak_window["start"] = 12:00:00 and peak_window["end"] = 21:00:00 (peak window is 9 hours, meaning that peak_window>n_control_window).

Should there be a minimum value on n_control_window? (not necessary - just asking!)

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Just checking - can peak_window['end'] be before peak_window['start']? Like if I wanted to define my peak window from 11pm (23:00:00) to 3am (03:00:00) - would that work?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No, there is no constraint on peak_window wrt n_control_window. n_control_window is just iterating over the 8760 is small chunks.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

_compute_peak_window_mask raises a ValueError if end < start.

performance_incentive (float)
Incentive revenue ($/kW per dispatch hour).
n_max_events (int)
Maximum discharge events per control window. Default 10.
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

great doc-strings! These are so helpful!

peak_window (dict)
Hours eligible for dispatch. 'start' and 'end' as HH:MM:SS strings.
performance_incentive (float)
Incentive revenue ($/kW per dispatch hour).
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
Incentive revenue ($/kW per dispatch hour).
Incentive revenue ($/commodity_rate_units per dispatch hour).

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I've only thought about this PR as a battery use case. We should chat about what this would mean.

m = pyomo.ConcreteModel(name="plm_dr")

# Parameters
P_max = self.config.max_charge_rate
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In order for the controller to be compatible with a case where we want to optimize the storage charge rate or capacity, then any method called in compute() should use inputs["max_charge_rate"] instead of self.config.max_charge_rate and inputs["storage_capacity"] instead of self.config.max_capacity.

An example of how to test this is somewhat shown in h2integrate/storage/test/test_storage_performance_model.py::test_generic_storage_with_simple_control_charge_rate_lessthan_demand, using the prob.set_val function (happy to chat about it if you want)

@vijay092 vijay092 force-pushed the peakload-optimized branch from b82857b to 410e488 Compare April 22, 2026 21:44
@johnjasa johnjasa self-requested a review April 23, 2026 19:45
Copy link
Copy Markdown
Collaborator

@jaredthomas68 jaredthomas68 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is great Sanjana! I'm excited about where this is going. I gave a lot of comments, but some of them may be irrelevant due to my experience with pyomo. The doc page details were very helpful, thank you.

**Given:**
- $\lambda_t$ := `supervisory_signal`: price, demand, or price $\times$ demand time series at time $t$
- $\mathcal{W}$ := `peak_window`: set of hours eligible for dispatch (e.g., 12:00--19:00)
- $\gamma$ := performance incentive (\$/kW per dispatch hour)
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is this truly required to be per hour, or is it per time-step? In other words, could it be per-minute if the user requested a smaller time step?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I was under the impression that H2I only supports hourly timesteps. I can make this more general.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We are currently in the process of enabling smaller time steps, so we are building the capability into new additions to the code!

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh, that's awesome!

- $\overline{\text{SoC}}$ := `max_soc_fraction`, $\quad \underline{\text{SoC}}$ := `min_soc_fraction`
- `n_control_window` := Horizon length for optimization
- $\mathcal{T} := \{0, 1, \ldots, T\}$: hourly time steps over `n_control_window`
- $\mathcal{M}_m$ := set of hours in month $m$, for $m = 1, \ldots, 12$
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same as above, does this have to be hours?

- $\eta_c$ := `charge_efficiency`, $\quad \eta_d$ := `discharge_efficiency`
- $\overline{\text{SoC}}$ := `max_soc_fraction`, $\quad \underline{\text{SoC}}$ := `min_soc_fraction`
- `n_control_window` := Horizon length for optimization
- $\mathcal{T} := \{0, 1, \ldots, T\}$: hourly time steps over `n_control_window`
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hours or time steps?


- $u_t \in \{0, 1\}$ := discharge binary: 1 if battery dispatches at hour $t$, 0 otherwise
- $v_t \in \{0, 1\}$ := charge binary: 1 if battery charges at hour $t$, 0 otherwise

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

hours or time steps?

Maximize total annual incentive revenue:

$$
\max_{u_t,\, v_t} \quad \gamma \cdot \bar{P} \sum_{t \in \mathcal{T}} u_t
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

should the objective be based on max charge/discharge rate, or on actual discharge rate at each time step?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I was assuming that the battery always charges and discharges at max discharge rate. If we want to solve for the actual discharge rate as well as when (the binary var) we want to discharge, the math becomes a little harder. I can work on that tomorrow.

)

m.objective = pyomo.Objective(
expr=-incentive * P_max * sum(m.discharge[t] for t in m.T),
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is there a good way to use possible discharge instead of rating (see my similar comment in docs)

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Totally fair point, I am working on it.

Comment on lines +418 to +419
+ eta_c * mdl.charge[t] * P_max / E_max
- mdl.discharge[t] * P_max / (eta_d * E_max)
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

see my comments in docs

RuntimeError: If GLPK returns a non-OK status or an
unacceptable termination condition.
"""
from pyomo.opt import SolverStatus, TerminationCondition
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

why is the import statement inside the method? Can you move it to the header?

pyomo.value(self.dr_model.objective),
)

def compute(self, inputs, outputs, discrete_inputs, discrete_outputs):
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it generally makes sense to keep the standard methods (init, setup, and compute) at the top of the class definition and have all other methods come after.

list[float]: ``(u_t - v_t) * P_max`` for each timestep in
the solved window. Positive = discharge, negative = charge.
"""
P_max = self.config.max_charge_rate
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think P_max should be adjusted to not overcharge the battery, but still allow it to go to max charge. Maybe you have the handled somehow and I missed it.

@vijay092 vijay092 force-pushed the peakload-optimized branch from 3d4735f to 2f0a325 Compare April 26, 2026 02:43
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants